/** https://ideone.com/VCu4IW **/
/*

Last Update date : 2020-10-30

直线交点
判断点在线段上
判断线段相交
求 p1->q 在 p1 -> p2 上投影的向量，即垂足坐标
镜像
点到线段
线段到线段最小值
p1 与 p2 之间弧度

多边形 有向 面积
射线法判点在多边形内 2:inside 1:on_seg, 0:out_side
凸包
凸包直径
q1->q2 左侧切出的凸包
ps 逆时针, 逆时针 a -> b, 求两个切点
凸包内切线(未实现?)
凸包外切线(未实现?)
最近点对

两圆之间关系
求圆与直线交点, 相切返回2个
求圆和圆的交点
求圆与圆外一点的切线
外切线
内切线
o圆心在原点，三角形1顶点位于原点的有向面积交
三角形外接圆圆心
三角形内接圆圆心, 配合 proj 求半径
最小包围圆随机 O(n)

半平面交
极角比较
凸包之间的最小距离
角度的并
角度的交

多边形内部最长距离
简单多边形三角剖分
快速检查线段是否和多边形严格相交
射线:st -> st+dir 与简单多边形的第一个严格交点
多边形 A 和直线 ( 线段 )k1->k2 严格相交
*/

#include <bits/stdc++.h>
#define rep(i,a,b) for(int i=(a);i<=(b);++i)
#define per(i,a,b) for(int i=(a);i>=(b);--i)
#define pb push_back
using namespace std;

typedef double db;

const db EPS = 1e-8;

inline int sgn(db a, db E=EPS) { return a < -E ? -1 : a > E; }
inline int cmp(db a, db b, db E=EPS) { return sgn(a-b,E); }

struct P {
  db x, y;
  P(){}
  P(db _x, db _y) : x(_x), y(_y) {}
  bool operator<(P o) const { return cmp(x,o.x) != 0? cmp(x,o.x) == -1 : cmp(y,o.y) == -1; }
  bool operator>(P o) const { return cmp(x,o.x) != 0? cmp(x,o.x) == 1 : cmp(y,o.y) == 1; }
  bool operator==(P o) { return cmp(x,o.x) == 0 && cmp(y,o.y) == 0; }
  P operator+(P p) { return P(x + p.x, y + p.y); }
  P operator-(P p) { return P(x - p.x, y - p.y); }
  P operator-() { return {-x,-y}; }
  P operator*(db d) { return P(x * d, y * d); }
  P operator/(db d) { return P(x / d, y / d); }
  db dot(P p) { return x * p.x + y * p.y; }
  db det(P p) { return x * p.y - y * p.x; }
  db disTo(P p) { return (*this - p).abs(); }
  db alpha() { return atan2(y, x); }
  void read() { cin >> x >> y; } // 根据题目修改
  db abs() { return sqrt(abs2()); }
  db abs2() { return x * x + y * y; }
  P rot90() { return P(-y, x); } // 逆时针
  P rot(db a) { db c = cos(a), s = sin(a); return {c*x-s*y, s*x+c*y}; } // 绕原点 a 弧度
  P unit() { return *this/abs(); }
  int quad() const { return sgn(y) == 1 || (sgn(y) == 0 && sgn(x) >= 0); } // 判断是否在上半区
};

struct L {
  P ps[2]; // ps[0] -> ps[1]
  L(){}
  L(P p1,P p2){ps[0]=p1, ps[1] = p2;}
  P& operator[](int i) { return ps[i]; }
  P dir() { return ps[1] - ps[0]; }
  db alpha() { return dir().alpha(); }
  bool include(P p) { return sgn((ps[1]-ps[0]).det(p - ps[0])) > 0; }
};

#define cross(p1,p2,p3) ((p2.x-p1.x)*(p3.y-p1.y)-(p3.x-p1.x)*(p2.y-p1.y))
#define crossOp(p1,p2,p3) sgn(cross(p1,p2,p3))

P isLL(P p1, P p2, P q1, P q2) { // 直线交点
  db a1 = cross(q1, q2, p1), a2 = -cross(q1, q2, p2);
  return (p1 * a2 + p2 * a1) / (a1 + a2);
}

P isLL(L l1,L l2){ return isLL(l1[0],l1[1],l2[0],l2[1]); }

bool isMiddle(db a, db m, db b) {
  return sgn(a - m) == 0 || sgn(b - m) == 0 || (a < m != b < m);
}

bool isMiddle(P a, P m, P b) { // 点在a b 之间的矩形内
  return isMiddle(a.x, m.x, b.x) && isMiddle(a.y, m.y, b.y);
}

bool onSeg(P p1, P p2, P q) { // 判断点在线段上
  return crossOp(p1,p2,q) == 0 && isMiddle(p1, q, p2);
}

bool onSeg_strict(P p1, P p2, P q) { // 判断点在线段上
  return crossOp(p1,p2,q) == 0 && sgn((q-p1).dot(p1-p2)) * sgn((q-p2).dot(p1-p2)) < 0;
}

bool intersect(db l1, db r1, db l2, db r2) {
  if(l1 > r1) swap(l1, r1); if(l2 > r2) swap(l2, r2);
  return !( cmp(r1,l2) == -1 || cmp(r2,l1) == -1 );
}

bool isSS(P p1, P p2, P q1, P q2) { // 判断线段相交，包含端点
  return intersect(p1.x,p2.x,q1.x,q2.x) && intersect(p1.y,p2.y,q1.y,q2.y) &&
  crossOp(p1,p2,q1) * crossOp(p1,p2,q2) <= 0 && crossOp(q1,q2,p1) * crossOp(q1,q2,p2) <= 0;
}

bool isSS_strict(P p1, P p2, P q1, P q2) { // 判断线段相交，不 包含端点
  return crossOp(p1,p2,q1) * crossOp(p1,p2,q2) < 0 && crossOp(q1,q2,p1) * crossOp(q1,q2,p2) < 0;
}

P proj(P p1, P p2, P q) { // 求 p1->q 在 p1 -> p2 上投影的向量，即垂足坐标，整点坐标注意，有可能垂足是小数
  P dir = p2 - p1;
  return p1 + dir * (dir.dot(q-p1) / dir.abs2());
}

P reflect(P p1, P p2, P q) { //镜像
  return proj(p1,p2,q)*2 - q;
}

db nearest(P p1, P p2, P q) { //点到线段
  P h = proj(p1, p2, q);
  if(isMiddle(p1,h,p2)) return h.disTo(q);
  return min(p1.disTo(q), p2.disTo(q));
}

db disSS(P p1, P p2, P q1, P q2) { // 线段到线段最小值
  if(isSS(p1,p2,q1,q2)) return 0;
  return min( min(nearest(p1,p2,q1), nearest(p1,p2,q2)), min(nearest(q1,q2,p1), nearest(q1,q2,p2)) );
}

db rad(P p1,P p2) { // p1 与 p2 之间弧度
  return atan2l(p1.det(p2),p1.dot(p2));
}

// polygon

db area(vector<P> ps) { // 多边形 有向 面积
  db ret = 0;
  rep(i,0,(int)ps.size()-1) ret+=ps[i].det(ps[(i+1)%ps.size()]);
  return ret/2;
}

int contain(vector<P> &ps, P p) { // 射线法判点在多边形内 2:inside 1:on_seg, 0:out_side
  int n = ps.size(), ret = 0;
  rep(i,0,n-1) {
  	P u = ps[i], v = ps[(i+1)%n];
  	if(onSeg(u, v, p)) return 1;
  	if(cmp(u.y, v.y) <= 0) swap(u, v);
  	if(cmp(p.y,u.y) > 0|| cmp(p.y,v.y) <= 0) continue;
  	ret ^= (crossOp(p,u,v) > 0);
  }
  return ret * 2;
}

vector<P> convexHull(vector<P> ps) { // 相邻边夹角小于180
  int n = ps.size(); if(n <= 1) return ps;
  sort(ps.begin(), ps.end());
  vector<P> qs(n * 2); int k = 0;
  for (int i = 0; i < n; qs[k++] = ps[i++])
  	while (k > 1 && crossOp(qs[k - 2], qs[k - 1], ps[i]) <= 0) --k;
  for (int i = n - 2, t = k; i >= 0; qs[k++] = ps[i--])
    while (k > t && crossOp(qs[k - 2], qs[k - 1], ps[i]) <= 0) --k;
  qs.resize(k - 1);
  return qs;
}

vector<P> convexHullNonStrict(vector<P> ps) { // 相邻边夹角小于等于180
  //需要预先对ps去重
  int n = ps.size(); if(n <= 1) return ps;
  sort(ps.begin(), ps.end());
  vector<P> qs(n * 2); int k = 0;
  for (int i = 0; i < n; qs[k++] = ps[i++])
		while (k > 1 && crossOp(qs[k - 2], qs[k - 1], ps[i]) < 0) --k;
  for (int i = n - 2, t = k; i >= 0; qs[k++] = ps[i--])
  	while (k > t && crossOp(qs[k - 2], qs[k - 1], ps[i]) < 0) --k;
  qs.resize(k - 1);
  return qs;
}

db convexDiameter(vector<P> ps) { // 旋转卡壳 凸包直径
  int n = ps.size(); if(n <= 1) return 0;
  int is = 0, js = 0, i, j;
  rep(i,0,n-1) is = ps[is]<ps[i]?is:i, js = ps[js]<ps[i]?i:js;
  db ret = ps[is].disTo(ps[js]);
  i = is, j = js;
  do {
    if( (ps[(i+1)%n]-ps[i]).det(ps[(j+1)%n]-ps[j]) >= 0) (++j)%=n;
    else (++i)%=n;
    ret = max(ret, ps[i].disTo(ps[j]));
  } while(i!=is || j!=js);
  return ret;
}

vector<P> convexCut(vector<P> ps, P q1, P q2) { // q1->q2 左侧切出的凸包
  int n = ps.size(), k = 0; if(n <= 1) return ps;
  vector<P> qs;
  rep(i,0,n-1) {
    P p1 = ps[i], p2 = ps[(i+1)%n];
    int c1 = crossOp(q1, q2, p1), c2 = crossOp(q1, q2, p2);
    if(c1 >= 0) qs.pb(p1); // left or onSeg
    if(c1*c2 < 0) qs.pb(isLL(q1,q2,p1,p2));
  }
  return qs;
}

void get_tangent(vector<P> ps, P p, int &a, int &b) { // ps 逆时针, 逆时针 a -> b, 求两个切点
  int n = ps.size();
  rep(i, 0, n-1) {
    if( crossOp(p,ps[i],ps[(i+1)%n]) >= 0 && crossOp(p, ps[(i+1)%n], ps[(i+2)%n]) < 0 ) a = (i+1)%n;
    if( crossOp(p,ps[i],ps[(i+1)%n]) < 0 && crossOp(p, ps[(i+1)%n], ps[(i+2)%n]) >= 0 ) b = (i+1)%n;
  }
}

db min_dist(vector<P>&ps,int l,int r){ // sort first 最近点对
  if(r-l+1<=5){
    db ret = 1e100;
    rep(i,l,r) rep(j,l,i-1) ret = min(ret,ps[i].disTo(ps[j]));
    return ret;
  }
  int m = (l+r)>>1;
  db ret = min(min_dist(ps,l,m),min_dist(ps,m+1,r));
  vector<P> qs; rep(i,l,r) if(abs(ps[i].x-ps[m].x)<= ret) qs.pb(ps[i]);
  sort(qs.begin(), qs.end(),[](P a,P b) -> bool {return a.y<b.y;}); /**c++11**/
  rep(i,1,(int)qs.size()-1) for(int j=i-1;j>=0&&qs[j].y>=qs[i].y-ret;--j) ret = min(ret,qs[i].disTo(qs[j]));
  return ret;
}

// circle

int type(P o1, db r1, P o2, db r2) { // 两圆之间关系
  db d = o1.disTo(o2);
  if(cmp(d, r1+r2)==1) return 4; // 相离
  if(cmp(d, r1+r2)==0) return 3; // 外切
  if(cmp(d, abs(r1-r2))==1) return 2; //相交
  if(cmp(d, abs(r1-r2))==0) return 1; // 内切
  return 0;
}

vector<P> isCL(P o,db r,P p1,P p2) { // 求圆与直线交点, 相切返回2个
  db x = (p1-o).dot(p2-p1), y = (p2-p1).abs2(), d = x * x - y * ((p1-o).abs2() - r*r);
  if(sgn(d) < 0) return {};
  d = max(d,0.0); P m = p1 - (p2-p1)*(x/y), dr = (p2-p1)*(sqrt(d)/y);
  return {m-dr,m+dr}; //along dir: p1->p2
}

vector<P> isCC(P o1, db r1, P o2, db r2) { // 需要检查是否是同一个圆，求圆和圆的交点，沿圆 o1 逆时针给出 , 相切给出两个
  db d = o1.disTo(o2);
  if(cmp(d, r1+r2) == 1) return {};
  d = min(d, r1+r2);
  db y = (r1 * r1 + d * d - r2 * r2) / (2 * d);
  db x = sqrt(r1 * r1 - y * y);
  P dr = (o2 - o1).unit();
  P q1 = o1 + dr * y, q2 = dr.rot90() * x;
  return {q1-q2,q1+q2};
}

vector<P> tanCP(P o, db r, P p) { // 求圆与圆外一点的切线, 沿圆 o 逆时针给出
  db x = (p - o).abs2(), d = x - r * r;
  if(sgn(d) <= 0) return {}; // 在圆上等于没切点
  P q1 = o + (p - o) * (r * r / x);
  P q2 = (p - o).rot90() * (r * sqrt(d) / x);
  return {q1-q2, q1+q2}; // 相对圆心逆时针
}

vector<L> extanCC(P o1, db r1, P o2, db r2) { // 外切线, 相对 o1 逆时针
  vector<L> ret;
  if(cmp(r1,r2) == 0) {
    P dr = (o2 - o1).unit().rot90() * r1;
    ret.pb({o1 + dr, o2 + dr}), ret.pb({o1 - dr, o2 - dr});
  }
  else {
    P p = (o2 * r1 - o1 * r2) / (r1 - r2);
    vector<P> ps = tanCP(o1,r1,p), qs = tanCP(o2,r2,p);
    rep(i,0,(int)min(ps.size(),qs.size())-1) ret.pb({ps[i],qs[i]});
  }
  return ret;
}

vector<L> intanCC(P o1, db r1, P o2, db r2) { // 内切线, 相对 o1 逆时针
  vector<L> ret;
  P p = (o1 * r2 + o2 * r1) / (r1 + r2);
  vector<P> ps = tanCP(o1,r1,p), qs = tanCP(o2,r2,p);
  rep(i,0,(int)min(ps.size(),qs.size())-1) ret.pb({ps[i],qs[i]});
  return ret;
}

db areaCT(db r, P p1, P p2) { // o圆心在原点，三角形1顶点位于原点的有向面积交
  vector<P> is = isCL(P(0,0), r, p1, p2); // 交点方向 p1 -> p2
  if(is.empty()) return r * r * rad(p1, p2) / 2;
  bool b1 = cmp(p1.abs2(),r*r) == 1, b2 = cmp(p2.abs2(),r*r) == 1;
  if(b1 && b2) {
    if(sgn((p1-is[0]).dot(p2-is[0])) <= 0)
      return r*r*(rad(p1,is[0]) + rad(is[1],p2))/2 + is[0].det(is[1])/2;
    else return r*r*rad(p1,p2)/2;
  }
  if(b1) return (r*r*rad(p1,is[0]) + is[0].det(p2))/2;
  if(b2) return (p1.det(is[1]) + r*r*rad(is[1],p2))/2;
  return p1.det(p2)/2;
}

P getOutCircle(P p1,P p2,P p3){ // 三角形外接圆圆心
  db a1=p2.x-p1.x,b1=p2.y-p1.y,c1=(a1*a1+b1*b1)/2;
  db a2=p3.x-p1.x,b2=p3.y-p1.y,c2=(a2*a2+b2*b2)/2;
  db d=a1*b2-a2*b1;
  return {p1.x+(c1*b2-c2*b1)/d,p1.y+(a1*c2-a2*c1)/d};
}

P getInCircle(P p1, P p2, P p3) { // 三角形内接圆圆心, 配合 proj 求半径
  db a=p2.disTo(p3), b=p1.disTo(p3), c=p2.disTo(p1);
  return {(a*p1.x+b*p2.x+c*p3.x)/(a+b+c),(a*p1.y+b*p2.y+c*p3.y)/(a+b+c)};
}

pair<P,db> min_circle(vector<P> ps){ // 最小包围圆随机 O(n)
  random_shuffle(ps.begin(),ps.end());
  int n = ps.size();
  P o = ps[0]; db r = 0;
  rep(i,1,n-1) if (o.disTo(ps[i]) > r + EPS) {
    o = ps[i], r = 0;
    rep(j,0,i-1) if (o.disTo(ps[j]) > r + EPS) {
      o = (ps[i]+ps[j])/2, r=o.disTo(ps[i]);
      rep(k,0,j-1) if (o.disTo(ps[k]) > r + EPS) {
        o = getOutCircle(ps[i],ps[j],ps[k]), r = o.disTo(ps[i]);
      }
    }
  }
  return {o,r};
}


// 半平面交
bool isInfPlane(vector<L>&l) { // 通过弧度覆盖的角度交，判断平面是否无穷
  vector<pair<db,db>> V;
  db pa, pb;
  for(int i = 0; i < l.size(); ++ i) {
    P dl = l[i].dir(), dr = -dl;
    db a = dl.alpha();
    db b = dr.alpha();
    if(a < 0) a += 2 * PI;
    if(b < 0) b += 2 * PI;
    if(i != 0) {
      pair<db,db> tt = {a,b}, ttt = {pa,pb};
      if( !angleIns(ttt,tt) ) return false; 
      pa = ttt.first, pb = ttt.second;
    }
    else pa = a, pb = b;
  }
  return true;
}

/** 某些情况下，当l0与l1交点足够远时，认为平行，以避免某些无用的交点带来的误差 **/
bool parallel(L l0, L l1) { 
  return sgn( l0.dir().det( l1.dir() ) ) == 0;
  // 交点坐标超出可能的范围时，认定为平行[-1e9,1e9]
  // P u = isLL(l0,l1);
  // return sgn( l0.dir().det( l1.dir() ) ) == 0 || cmp(abs(u.x),1e9) == 1 || cmp(abs(u.y),1e9) == 1; 
}
bool sameDir(L l0, L l1) { return parallel(l0, l1) && sgn(l0.dir().dot(l1.dir()) ) == 1; }
bool cmpa (P a,  P b) { // 极角比较
  if (a.quad() != b.quad()) { return a.quad() < b.quad(); }
  else { return sgn( a.det(b) ) > 0; }
}
bool operator < (L l0, L l1) {
  if (sameDir(l0, l1)) { return l1.include(l0[0]); }
  else { return cmpa( l0.dir(), l1.dir() ); }
}
// bool operator < (L l0, L l1) { /** 当数值范围较大(>1e9),精度较低,使用弧度排序,最好预处理弧度 **/
// 	if ( cmp(l0.alpha() , l1.alpha()) != 0 ) return cmpa(l0.dir(),l1.dir());
// 	return l1.include(l0[0]);
// }
bool check(L u, L v, L w) { return w.include(isLL(u,v)); }
vector<P> halfPlaneIS(vector<L> &l) { // 求半平面交 , 半平面是逆时针方向 , 输出按照逆时针,  注意数据范围一般会引入的4个边界
  sort(l.begin(), l.end());
  // if(isInfPlane(l)) {/*INF*/}
  deque<L> q;
  for (int i = 0; i < (int)l.size(); ++i) {
    if (i && sameDir(l[i], l[i - 1])) continue;
    // if (i && cmp(l[i].alpha(),l[i-1].alpha()) == 0) continue; /**使用弧度排序时**/
    while (q.size() > 1 && !check(q[q.size() - 2], q[q.size() - 1], l[i])) q.pop_back();
    while (q.size() > 1 && !check(q[1], q[0], l[i])) q.pop_front();
    q.push_back(l[i]);
  }
  while (q.size() > 2 && !check(q[q.size() - 2], q[q.size() - 1], q[0])) q.pop_back();
  while (q.size() > 2 && !check(q[1], q[0], q[q.size() - 1])) q.pop_front();
  vector<P> ret;
  for (int i = 0; i < (int)q.size(); ++i) ret.push_back(isLL(q[i], q[(i + 1) % q.size()]));
  return ret;
}



// O(nlogn) 半平面交求凸包的交面积
db getConvConv(vector<P> ps, vector<P> qs) {
  if(area(ps) < 0) reverse(ps.begin(), ps.end());
  if(area(qs) < 0) reverse(qs.begin(), qs.end());
  vector<L> Ls;
  int n = ps.size(), m = qs.size();
  rep(i, 0, n-1) Ls.pb(L(ps[i],ps[(i+1)%n]));
  rep(i, 0, m-1) Ls.pb(L(qs[i],qs[(i+1)%m]));
  vector<P> ret = halfPlaneIS(Ls);
  return area(ret);
}
// O(nlogn) 半平面交判凸包是否有交
int isConvConv(vector<P> ps, vector<P> qs) {	
  if( sgn(getConvConv(ps,qs)) ) return 1;
  return 0;
}
// O(n) 凸包之间的最小距离
// 凸包会相交时，加入注释部分
#define Pii pair<int,int>
#define fi first
#define se second
db disConvConv(vector<P> ps, vector<P> qs) {
  // if( isConvConv(ps,qs) ) return 0;
  int n = ps.size(), m = qs.size(), i = 0, j = 0, f = 0;
  vector<Pii> d = {{0,1},{0,m-1},{1,0},{n-1,0}};
  auto cal_d = [&](int u, int v) {
    return disSS(ps[u],ps[(u+1)%n], qs[v],qs[(v+1)%m]);
  };
  db ans = cal_d(i,j);
  while(true) { f = 0;
    rep(k, 0, 3) if(ans > cal_d((i+d[k].fi)%n,(j+d[k].se)%m)) {
      f = 1;
      while(ans > cal_d((i+d[k].fi)%n,(j+d[k].se)%m)) 
        i=(i+d[k].fi)%n, j=(j+d[k].se)%m, ans = cal_d(i,j);
    }
    if(!f) break;
  }
  return ans;
}

//角度的并
vector<pair<db,db>> angleUnion(vector<pair<db,db>> segs){ // [0,2*PI]
  if(segs.size() <= 1) return segs;
  sort(segs.begin(), segs.end());
  pair<db,db> seg = segs[0];
  vector<pair<db,db> > ans;
  for(int i = 1; i < segs.size(); ++ i) {
    auto s = segs[i];
    if(s.first < seg.second) seg.second = max(seg.second,s.second);
    else ans.push_back(seg), seg = s;
  }
  ans.push_back(seg);
  return ans;
}
vector<pair<db,db>> angleUnion(vector<pair<P,P>> segs){ // [pl,pr] pl->pr 逆时针
  vector<pair<db,db>> ans;
  for(int i = 0; i < segs.size(); ++ i) {
    db al = segs[i].first.alpha();
    db ar = segs[i].second.alpha();
    if(al < 0) al += 2*PI;
    if(ar < 0) ar += 2*PI;
    if(al < ar) ans.pb({al,ar});
    else ans.pb({al, 2*PI}), ans.pb({0, ar});
  }
  return angleUnion(ans);
}
//角度的交(待测)
// ta,tb 可以为锐角、钝角及平角，值域为[0,2*PI]，相交返回true，需要注意精度问题
bool angleIns(pair<db,db> &ta, pair<db,db> tb) {
  auto a = ta;
  auto b = tb;
  if(a.first > a.second) a.second += 2*PI;
  if(b.first > b.second) b.second += 2*PI;
  vector<pair<db,db>> sa, sb,ans;
  sa.emplace_back(a.first,a.second);
  sb.emplace_back(b.first,b.second);
  sa.emplace_back(a.first+2*PI,a.second+2*PI);
  sb.emplace_back(b.first+2*PI,b.second+2*PI);
  for(auto &A: sa)for(auto &B: sb) {
    db l = max(A.first,B.first);
    db r = min(A.second,B.second);
    if(l <= r) { // 角度eps比较敏感，gym/102428D测试中不设eps可通过
      if(l > 2*PI) l -= 2*PI;
      if(r > 2*PI) r -= 2*PI;
      ans.emplace_back(l,r);
    }
  }
  if(ans.size()>0) ta = ans[0];
  return (!ans.empty());
}


//简单多边形
db calLongestDisInPo(vector<P> ps) { // 多边形内部最长距离 含边界
  db ans = 0; int n = ps.size();
  function<void(P,P)> update = [&] (P p1, P p2) -> void {
    vector<P> qs;
    rep(i,0,n-1) if( sgn((ps[i]-ps[(i+1)%n]).det(p1-p2)) != 0 ) {
      P is = isLL(ps[i],ps[(i+1)%n],p1,p2); if( isMiddle(ps[i],is,ps[(i+1)%n]) ) qs.pb(is);
    }
    sort(qs.begin(),qs.end()); qs.erase(unique(qs.begin(),qs.end()),qs.end());
    if( cmp(qs[0].disTo(qs[qs.size()-1]), ans) <= 0 ) return;
    int cc = qs.size();
    int pos1 = lower_bound(qs.begin(),qs.end(),p1) - qs.begin();
    int pos2 = lower_bound(qs.begin(),qs.end(),p2) - qs.begin();
    if(pos1 > pos2) swap(pos1, pos2);
    rep(i,pos1+1,pos2) if(!contain(ps, (qs[i] + qs[i-1]) / 2 )) return;
    db tmp = p1.disTo(p2);
    per(i,pos1,1) { if(!contain(ps, (qs[i] + qs[i-1]) / 2)) break; tmp += qs[i].disTo(qs[i-1]); }
    rep(i,pos2+1,cc-1) { if(!contain(ps, (qs[i] + qs[i-1]) / 2)) break; tmp += qs[i].disTo(qs[i-1]); }
    ans = max(ans, tmp);
    return;
  };
  rep(i,0,n-1) rep(j,i+1,n-1) update(ps[i], ps[j]);
  return ans;
}



// O(n^2) 简单多边形三角剖分 多边形逆时针
typedef vector<vector<P>> VVP;
template <class A, class B> A pr(A it, B &L) { return it == L.begin() ? --L.end() : --it; };
template <class A, class B> A nx(A it, B &L) { return ++it == L.end() ? L.begin() : it; };
VVP divPtoTFast(vector<P> ps) {
#define IT list<P>::iterator
#define ITE list<IT>::iterator
  if( area(ps) < 0 ) reverse(ps.begin(), ps.end()); //保证逆时针
  list<P> L(ps.begin(), ps.end()); list<IT> E; int n = ps.size();
  auto isEar = [&](IT it) { // 注意crossOp 中 *it 加括号
    P B = *it, A = *pr(it, L), C = *nx(it, L);
    if (crossOp(B, C, A) < 0) return 0;
    for (IT it = L.begin(); it != L.end(); ++it) {
      IT u = it, v = nx(u, L);
      if (crossOp(A, B, (*u)) > 0 && crossOp(B, C, (*u)) > 0 && crossOp(C, A, (*u)) > 0) return 0;
      if (isSS_strict(*u, *v, A, C)) return 0;
    }
    return 1;
  };
  for (IT it = L.begin(); it != L.end(); ++it) if (isEar(it)) E.pb(it);
  ITE ite = E.begin(), ep, en; IT now, pp, nn;
  VVP ans(n - 2);
  for (int f, i = 0; i < n - 2; ++i) {
    now = *ite, pp = pr(now, L), nn = nx(now, L), ep = pr(ite, E), en = nx(ite, E), ans[i] = {*now, *pp, *nn}, L.erase(now);
    auto zypnb = [&](IT pp, ITE ep, ITE S) {
      if ((f = isEar(pp)) && *ep != pp) E.insert(S, pp);
      else if (!f && *ep == pp) E.erase(ep);
      return;
    }; zypnb(pp, ep, ite), zypnb(nn, en, en);
    if ((ite = E.erase(ite)) == E.end()) ite = E.begin();
  }
  return ans;
#undef ITE
#undef IT
}

/*  传说 O(n^3) 但是实际好像 O(n^2) 的 简单多边形三角剖分 多边形逆时针
template<class A,class B> A nx(A it,B &L){ return ++it == L.end() ? L.begin() : it; };
typedef vector<P> VP;
typedef vector<vector<P>> VVP;
VVP divPtoT(vector<P> ps){
  VVP T; T.clear();
  if(ps.size()<3) return T;
  list<int> L; rep(i, 0, ps.size()-1) L.pb(i);
  auto ck = [&](P A, P B, P C) {
    if( crossOp(B,C,A) < 0 ) return 0;
    for(auto it = L.begin(); it != L.end(); ++ it) {
      auto u = it, v = nx(u,L);
      if(crossOp(A,B,ps[(*u)])>0&&crossOp(B,C,ps[(*u)])>0&&crossOp(C,A,ps[(*u)])>0) return 0;
      if(isSS_strict(ps[*u],ps[*v], A,C)) return 0;
    }
    return 1;
  };
  P A,B,C;
  while(L.size() > 3) {
    auto it = L.begin();
    A = ps[*it]; ++it; B = ps[*it]; --it; C = ps[*--L.end()];
    if(ck(A,B,C)) { L.erase(it); T.pb({A,B,C}); continue; }
    auto ed = --L.end();
    B = ps[*ed]; -- ed; A = ps[*ed]; ++ ed; C = ps[*L.begin()];
    if(ck(A,B,C)) { L.erase(ed); T.pb({A,B,C}); continue; }
    it = ++L.begin();
    for( ; it != ed; ++ it) {
      B = ps[*it]; --it; A = ps[*it]; ++it; ++it; C = ps[*it]; --it;
      if(ck(A,B,C)) { L.erase(it); T.pb({A,B,C}); break; }
    }
  }
  if(L.size() == 3) {
    vector<P> tmp; 
    for(auto it = L.begin(); it != L.end(); ++ it) tmp.pb(ps[*it]);
    T.pb(tmp);
  }
  return T;
} */

int checkinp(P r,P l,P m){
  if (cmpa(l,r)){return cmpa(l,m)&&cmpa(m,r);}
  return cmpa(l,m)||cmpa(m,r);
}
int checkPosFast(vector<P> &A, P k1, P k2){ // 快速检查线段是否和多边形严格相交
  if (contain(A,k1)==2||contain(A,k2)==2) return 1; if (k1==k2) return 0;
  A.push_back(A[0]); A.push_back(A[1]);
  for (int i=1;i+1<A.size();i++)
    if ( cmp((k1-A[i-1]).det(k2-A[i-1]), (k1-A[i]).det(k2-A[i])) != 0  ) {
      P now=isLL(A[i-1],A[i],k1,k2);
      if (isMiddle(A[i-1],now,A[i])==0||isMiddle(k1,now,k2)==0) continue;
      if (now==A[i]){
        if (A[i]==k2) continue;
        P pre=A[i-1],ne=A[i+1];
        if (checkinp(pre-now,ne-now,k2-now)) {
          A.pop_back();A.pop_back();return 1;
        }
      } else if (now==k1){
        if (k1==A[i-1]||k1==A[i]) continue;
        if (checkinp(A[i-1]-k1,A[i]-k1,k2-k1)) {
          A.pop_back();A.pop_back();return 1;
        }
      } else if (now==k2||now==A[i-1]) continue;
      else {
        A.pop_back();A.pop_back();return 1;
      }
    }
  A.pop_back();A.pop_back();
  return 0;
}

P getRayPloy(vector<P> &ps, P st, P dir) { // 射线:st -> st+dir 与简单多边形的第一个严格交点 
  P res = st + dir.unit() * 1e5;
  rep(i,0,(int)ps.size()-1) {
  	P u = ps[i], v = ps[(i+1)%ps.size()];
  	if(sgn((res-st).det(u-v)) == 0) continue;
  	P tt = isLL(res,st,u,v);
  	if(onSeg_strict(res,st,tt) && onSeg(u,v,tt)) {
  	  res = tt;
  	}
  }
  if( contain(ps,(st+res)/2) ) return st;
  return res;
}


// 待测 code ...
int checkPoS(vector<P>A,P k1,P k2) { // 多边形 A 和直线 ( 线段 )k1->k2 严格相交 , 注释部分为线段
  struct ins{
    P m,u,v;
    int operator < (const ins& k) const {return m<k.m;}
  }; vector<ins>B;
  //if (contain(A,k1)==2||contain(A,k2)==2) return 1;
  vector<P>poly=A; A.push_back(A[0]);
  for (int i=1;i<A.size();i++) if ( cmp((k1-A[i-1]).det(k2-A[i-1]), (k1-A[i]).det(k2-A[i])) != 0 ) {
    P m = isLL(A[i-1],A[i],k1,k2);
    if (isMiddle(A[i-1],m,A[i])/*&&isMiddle(k1,m,k2)*/) B.push_back((ins){m,A[i-1],A[i]});
  }
  if (B.size()==0) return 0; sort(B.begin(),B.end());
  int now=1; while (now<B.size()&&B[now].m==B[0].m) now++;
  if (now==B.size()) return 0;
  int flag=contain(poly,(B[0].m+B[now].m)/2);
  if (flag==2) return 1;
  P d=B[now].m-B[0].m;
  for (int i=now;i<B.size();i++){
    if (!(B[i].m==B[i-1].m)&&flag==2) return 1;
    int tag=sgn( (B[i].v-B[i].u).det(B[i].m+d-B[i].u) );
    if (B[i].m==B[i].u||B[i].m==B[i].v) flag+=tag; else flag+=tag*2;
  }
  //return 0;
  return flag==2;
}

// 圆并/交
// 三维

int main() {
  return 0;
}
